Part I: Data Description¶
The dataset consists of 17 variables with a sample space of 252 men. The 17 variables are INDO (index), Percentage of body fat (%), Body density from underwater weighing (gm/cm^3), Age (year), Weight (lbs), Height (inches), Adioposity (bmi) and ten Body Circumferences (Neck, Chest, Abdomen, Hip, Thigh, Knee, Ankle, Biceps, Forearm, Wrist, all in units of cm). Percentage of body fat is given from Siri's (1956) equation:
$$BodyFat\ \% = \frac{495}{Density}\ –\ 450$$
Part II: Thesis Statement¶
To accurately estimate bodyfat with clinical measurements, we used several criterions to select important variables. After model diagonsis, we transformed independent variables to better satisfy linear model assumptions.
We found out that among all 14 variables, the linear function of WEIGHT, WEIGHT transformation, ABDOMEN, FOREARM, and WRIST can best interpret bodyfat.
library("MASS")
library("car")
library(plotly)
data = read.csv("Bodyfat.csv", header = TRUE)[,-1]
attach(data)
plot_ly(data, x = ~BODYFAT, type="histogram")
We can see that there is one point that 0% bodyfat, which is impossible; also, the maximum value of bodyfat is abnormal. Now, let's check these two points:
data[which(BODYFAT == 0), ]
data[which(BODYFAT == max(data$BODYFAT)), ]
bad <- NULL
bad <- which(BODYFAT == 0)
We tried to predict bodyfat No.182 with density by Siri's equation, but the prediction is negative, so we considered the observation has missing dependent veriable, so we rule it out from model.
In addition, the Maximums of Weight, Neck Cir., Chest Cir., Abdomen Cir., Hip Cir., Thigh Cir., Knee Cir., Biceps Cir., Wrist Cir. are all from the same observation, which is the 39th. We now set it aside for further study.
We conclude that the record of 182 is not valid.
General Demographic Information¶
We fit full linear model without No.182 record:
m1 <- lm(BODYFAT ~ ., data = data[-bad, -2])
plot(m1, which = 4)
abline(h = 4/(nrow(data)-ncol(data)), lty = 2)
data[c(39, 42, 86),]
bad <- c(bad,42)
We conclude that the record of 42 may be not valid.
Consistence of DENSITY versus BODYFAT¶
reverse_de <- 1/DENSITY
plot_ly(data, x = ~reverse_de, y = ~BODYFAT, type= "scatter", mode = "markers")
m0 <- lm(BODYFAT ~ reverse_de)
plot(m0, which = 1)
data[96,]
round(495/DENSITY[96]- 450, 2)
data[c(76, 24),]
round(495/DENSITY[76]- 450, 2)
data[c(48, 24),]
round(495/DENSITY[48]- 450, 2)
data$BODYFAT[48] <- round(495/data$DENSITY[48]- 450, digits = 1)
We decided to use the BODYFAT value calculated with DENSITY.
Consistence of BMI versus WEIGHT & HEIGHT¶
BMI <- (WEIGHT/2.2046226218)/((HEIGHT*0.0254)^2)
boxplot(BMI - ADIPOSITY)
which(abs(BMI - data$ADIPOSITY) > 1)
data[which(WEIGHT > 183 & WEIGHT < 185 & HEIGHT > 67 & HEIGHT < 69), ]
data[which(WEIGHT > 153 & WEIGHT < 155 & HEIGHT > 69 & HEIGHT < 71), ]
We conclude that the record of 163, 221 may be not valid.
bad <- c(bad, 221, 163)
detach(data)
data <- data[-bad, -2]
attach(data)
m1 <- lm(BODYFAT ~ ., data)
m_null <- lm(BODYFAT ~ 1, data)
Variable Selection¶
| Method | Selected Varibles |
|---|---|
| BIC Backward | WEIGHT, ABDOMEN, FOREARM, WRIST |
| BIC Forward & Both | ABDOMEN, WEIGHT |
| AIC Backward | 10 variables |
| AICForward & Both | 6 variables |
| Mallow's Cp | 9 variables |
| LASSO | 5 variables |
### AIC
m_null <- lm(BODYFAT ~ 1, data)
m_AIC_back <- step(m1, k=2)
m_AIC_for <- step(m_null, direction="forward",
scope=list(lower=~1,upper=m1))
m_AIC_both <- step(m_null, direction="both",
scope=list(lower=~1, upper=m1)) # the selected models seem to complicated
m_BIC_back <- step(m1, k=log(nrow(data)-1)) # WEIGHT, ABDOMEN, FOREARM, WRIST
m_BIC_for <- step(m_null, direction="forward",
scope=list(lower=~1,upper=m1), k=log(nrow(data)-1)) # WEIGHT, ABDOMEN
m_BIC_both <- step(m_null, direction="both",
scope=list(lower=~1,upper= m1), k=log(nrow(data)-1)) # WEIGHT, ABDOMEN
m2 <- m_BIC_both # keep only ABDOMEN, WEIGHT
(s2 <- summary(m2))
(mse2 <- sum((s2$residuals)^2)/nrow(data))
# round the model to make it eaasier to calculate
fit <- -40 + ABDOMEN - 0.2*WEIGHT
mse <- sum((fit - BODYFAT)^2)/nrow(data)
res <- fit - BODYFAT
m3 <- lm(BODYFAT ~ ABDOMEN)
summary(m3)
m4 <- lm((BODYFAT)*WEIGHT ~ ABDOMEN + WEIGHT, data) # transform
summary(m4)
(mse4 <- sum((m4$residuals/WEIGHT)^2)/nrow(data)) # worse
m5 <- lm(BODYFAT ~ WEIGHT + ABDOMEN + FOREARM + WRIST, data) # WEIGHT, ABDOMEN, FOREARM, WRIST
(s5 <- summary(m5))
(mse5 <- sum((s5$residuals)^2)/nrow(data))
# round the model
fit5 <- -35 -0.15 * WEIGHT + ABDOMEN + 0.4*FOREARM - WRIST
mse5 <- sum((fit5- BODYFAT)^2)/nrow(data) # harder to round
anova(m5, m1)
We did ANOVA test on Full model versus BIC backward model (4 variables), the null hypothesis was retained.
anova(m2, m5)
We did ANOVA test on BIC backward model (4 variables) versus BIC forward model (2 variables), the null hypothesis was rejected.
We decided to reserve the two models mentioned above for further improvement.
Model Diagnosis¶
Two-variables Models¶
This model passed all of the normality test, multi-collinearity test and homoscedasticity test.
Test Linearity:
crPlots(m2)
bc <- boxcox(WEIGHT~1, data = data, lambda = seq(-10, 10, length = 10))
trans <- bc$x[which.max(bc$y)]
W2 <- WEIGHT^trans
mt<- lm(BODYFAT ~ ABDOMEN + WEIGHT + W2)
par(mfrow = c(2,2))
plot(mt)
par(mfrow = c(1,1))
crPlots(mt)
Four-variables Models¶
This model passed all of the normality test, multi-collinearity test and homoscedasticity test.
crPlots(m5)
mt_2 = lm(BODYFAT~WEIGHT+W2+ABDOMEN+FOREARM+WRIST, data)
Linearity Test
crPlots(mt_2)
Normality Test
shapiro.test(mt_2$residuals)
new <- data.frame(BODYFAT = BODYFAT, final = mt_2$fitted.value)
orders = order(new$BODYFAT)
plot_ly(new) %>%
add_trace(x = ~1:248, y = ~new$BODYFAT[orders], type = "scatter", color = "True",
marker = list(color = '#E69F00'), mode = "marker") %>%
add_trace(x = ~1:248, y = ~new$final[orders], type = "scatter", color = "Estimated",
marker = list(color = '#56B4E9'),mode = "marker") %>%
layout(title = 'True BodyFat Values v.s. Estimated',
xaxis = list(title = "Index in Increasing Ordered"),
yaxis = list(title = "BodyFat"))
Multi-collinearity Test
vif(mt_2)
Homoscedasticity Test
ncvTest(mt_2)
Proposed SLR Model¶
$$BodyFat\ \% = 50.183 - 0.2571Weight - 720.904Weight^{-0.505} + 0.904Abdomen + 0.279Forearm - 1.308Wrist$$
Strengths and Weaknesses of the Analysis¶
The model is a reasonable model between body fat % and abdomen, weight.
Moreover, the model has the following strengths and advantages:
- Linearity: seems reasonable. According the Component & Partial Residual plot and model residule v.s. fitted value plot.
- Explanatory variables: reasonable since it's easy to think one’s body fat with his weight and abdomen. Because intuitively, if a man has a big abdomen, he tends to be fatter. Meanwhile, if his weight is heavier given other body measurements fixed, he is likely to be more muscled, since fat has less density. Besides, the variable forearm and wrist can be treated as a measurement of body frame size which will also contribute to body fat's calculation
- Constant effects: reasonable because it is not affected by age or non-body factors.
- Normally distributed errors: seems reasonable from the Q-Q plot diagnostic and also the model passed Shapiro-Wilk’s test.
- Constant variance: the model passed the score test for non-constant error variance
Overall, our model provides a relatively simple way of predict the body fat % purely based on weight, forearm, wrist and abdomen.
Last but not least, there still exists some potential weaknesses or questions.
- Should there be non-linear relationships? We only fit linear models and the true relations could be more complicated.
- How to solve the multicollinearity between circumferences and other variables?
- Could the conclusion also be used for women? The data was only collected for men, and therefore the model is only suitable for men. Does there exist a general formula for both men and women?
Rule of Thumb¶
- "multiply your forearm circumference (cm) by 0.5, add the difference of your abdomen circumference (cm) and wrist circumference (cm), minus one length of your weight (lb) and minus 35"
- "multiply your weight (lb) by 0.2, add your abdomen (cm) and minus 40"
Example Usage:¶
for a 170lbs man with abdomen circumference about 90 cm, forearm circumference about 28 cm and wrist circumference about 17 cm, his predicted body fat % percentage would be around 19.55%. There is a 95% probability that his body fat is between 18.54% and 20.55%.
With the rule of thumb, you get about 16% as the predicted body fat %.